{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Ex3: Fitting spiky responses"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The Vector Fitting feature is demonstrated using a 4-port example network copied from the scikit-rf `tests` folder. This network is a bit tricky to fit because of its many resonances in the individual response. Additional explanations and background information can be found in the [Vector Fitting tutorial](../../tutorials/VectorFitting.ipynb)."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "import matplotlib.pyplot as mplt\n",
    "import numpy as np\n",
    "\n",
    "import skrf"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "To create a `VectorFitting` instance, a `Network` containing the frequency responses of the N-port is passed. In this example a copy of `Agilent_E5071B.s4p` from the `skrf/tests` folder is used:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "nw = skrf.network.Network('./Agilent_E5071B.s4p')\n",
    "vf = skrf.VectorFitting(nw)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For comparison, this network will be fitted using both the regular vector fitting routine `vector_fit()` as well as the automatic vector fitting routine `auto_fit()`."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Regular Vector Fitting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For regular vector fitting, the number and type of the initial poles have to be specified, which depend on the *behaviour* of the responses. As a rule of thumb for an initial guess, one can count the number of resonances or \"bumps\" in the individual responses. In this case, the 4-port network has 16 responses to be fitted. As shown in the magnitude plots below, $S_{11}$ and some other responses are quite *spiky* and have roughly 15 local maxima each and about the same number of local minima in between. Other responses have only 5-6 local maxima, or they are very noisy with very small magnitudes (like $S_{24}$ and $S_{42}$). Assuming that most of the 15 maxima of $S_{11}$ occur at the same frequencies as the maxima of the other responses, one can expect to require 15 complex-conjugate poles for a fit. As this is probably not completely the case, trying with 20-30 poles should be a good start to fit all of the resonances in all of the responses."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# plot magnitudes of all 16 responses in the 4-port network\n",
    "fig, ax = mplt.subplots(4, 4)\n",
    "fig.set_size_inches(12, 8)\n",
    "for i in range(4):\n",
    "    for j in range(4):\n",
    "        nw.plot_s_mag(i, j, ax=ax[i][j])\n",
    "        ax[i][j].get_legend().remove()\n",
    "fig.tight_layout()\n",
    "mplt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "After trying different numbers of real and complex-conjugate poles, the following setup was found to result in a very good fit. The initial model order is $1 + 2 * 26 = 53$, which will remain unchanged during the fitting process. Other setups were also found to work well for this network (e.g. 0-2 real poles and 25-26 cc poles):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vf.vector_fit(n_poles_real=1, n_poles_cmplx=26)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f'model order = {vf.get_model_order(vf.poles)}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The convergence can also be checked with the convergence plot:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vf.plot_convergence()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The fitted model parameters are now stored in the class attributes `poles`, `residues`, `proportional_coeff` and `constant_coeff` for further use. To verify the result, the fitted model responses can be compared to the original network responses. As the model will return a response at any given frequency, it makes sense to also check its response outside the frequency range of the original samples. In this case, the original network was measured from 0.5 GHz to 4.5 GHz, so we can plot the fit from dc to 10 GHz:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "freqs = np.linspace(0, 10e9, 501)\n",
    "fig, ax = mplt.subplots(4, 4)\n",
    "fig.set_size_inches(12, 8)\n",
    "for i in range(4):\n",
    "    for j in range(4):\n",
    "        vf.plot_s_mag(i, j, freqs=freqs, ax=ax[i][j])\n",
    "        ax[i][j].get_legend().remove()\n",
    "fig.tight_layout()\n",
    "mplt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f'rms error = {vf.get_rms_error()}')"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "As shown in the plots, a very good fit was achieved. This is also indicated with the low rms error of less than 0.01. However, a UserWarning about a non-passive fit was printed (see output of `vector_fit()` above). The assessment and enforcement of model passivity is described in more detail in [this example](./vectorfitting_ex4_passivity.ipynb), which is important for certain use cases in circuit simulators, i.e. transient simulations. To use the model in a circuit simulation, an equivalent circuit can be created based on the fitting parameters. This is currently only implemented for SPICE, but the structure of the equivalent circuit can be adopted to any kind of circuit simulator."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "`vf.write_spice_subcircuit_s('/home/vinc/Desktop/4-port_model.sp')`"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The exported `.sp` file can then be imported into SPICE as a subcircuit. Have a look at the [Ring Slot Example](./vectorfitting_ex1_ringslot.ipynb)."
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Automatic Vector Fitting"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For the following automatic vector fitting routine, the user does not have to specify the starting poles. The default settings should be sufficient for a successful fit, as the required number of poles gets optimized during the process. Still, changes in the parameters for `auto_fit()` can further improve the convergence and reduce the model order."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "vf.auto_fit()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f'model order = {vf.get_model_order(vf.poles)}')\n",
    "print(f'n_poles_real = {np.sum(vf.poles.imag == 0.0)}')\n",
    "print(f'n_poles_complex = {np.sum(vf.poles.imag > 0.0)}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "print(f'rms error = {vf.get_rms_error()}')"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "freqs = np.linspace(0, 10e9, 501)\n",
    "fig, ax = mplt.subplots(4, 4)\n",
    "fig.set_size_inches(12, 8)\n",
    "for i in range(4):\n",
    "    for j in range(4):\n",
    "        vf.plot_s_mag(i, j, freqs=freqs, ax=ax[i][j])\n",
    "        ax[i][j].get_legend().remove()\n",
    "fig.tight_layout()\n",
    "mplt.show()"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3 (ipykernel)",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.11.2"
  },
  "nbsphinx": {
   "timeout": 240
  }
 },
 "nbformat": 4,
 "nbformat_minor": 4
}
